mashR Alternatives

Goal: get more interesting patterns from mashR

library(flashr)

Error in library(flashr) : there is no package called ‘flashr’

Get QTL info

data.strong = mash_set_data(beta.hat[top.hits,], se.hat[top.hits,], V=Vhat, df=df.hat[top.hits,])

the condition has length > 1 and only the first element will be used

What do the existing mashR matrices look like?

names(m$fitted_g$Ulist)

 [1] "ED_PCA_1"                   "ED_PCA_2"                   "ED_PCA_3"                   "ED_PCA_4"                  
 [5] "ED_PCA_5"                   "ED_tPCA"                    "identity"                   "Astrocytes"                
 [9] "Visceral-neurons"           "Amacrine-cells"             "Oligodendrocytes"           "Ganglion-cells"            
[13] "Smooth-muscle-cells"        "Limbic-system-neurons"      "Vascular-endothelial-cells" "Photoreceptor-cells"       
[17] "CCL19-CCL21-positive-cells" "CLC-IL5RA-positive-cells"   "Mesangial-cells"            "AFP-ALB-positive-cells"    
[21] "Bipolar-cells"              "Satellite-cells"            "ENS-neurons"                "ENS-glia"                  
[25] "Granule-neurons"            "Hepatoblasts"               "Microglia"                  "IGFBP1-DKK1-positive-cells"
[29] "Purkinje-neurons"           "Epicardial-fat-cells"       "Cardiomyocytes"             "equal_effects"             
[33] "simple_het_1"               "simple_het_2"               "simple_het_3"              

Here we've got 5 "ED PCA" matrices, one "ED tPCA" matrix, the identity matrix, specific effects for each context, equal effects, and simple het effects. What do the ED PC matrices look like?

corrplot(m$fitted_g$Ulist$ED_PCA_1)

corrplot(m$fitted_g$Ulist$ED_PCA_2)

corrplot(m$fitted_g$Ulist$ED_PCA_3)

corrplot(m$fitted_g$Ulist$ED_PCA_4)

corrplot(m$fitted_g$Ulist$ED_PCA_5)

corrplot(m$fitted_g$Ulist$ED_tPCA)

Sooo that's a problem! How about those het ones?

corrplot(m$fitted_g$Ulist$simple_het_1)

corrplot(m$fitted_g$Ulist$simple_het_2)

corrplot(m$fitted_g$Ulist$simple_het_3)

So these... all suck, no wonder we are not seeing interesting patterns. Looking under the hood, lets see if maybe we're running into some issues just from PCA

U.pca = cov_pca(data.strong, 5)
corrplot(U.pca$PCA_1, is.corr=F)

corrplot(U.pca$PCA_2, is.corr=F)

corrplot(U.pca$PCA_3, is.corr=F)

corrplot(U.pca$PCA_4, is.corr=F)

corrplot(U.pca$PCA_5, is.corr=F)

corrplot(U.pca$tPCA, is.corr=F)

If we want to see the patterns coming out, we obviously need to have some patterns going in, and right now, we do not.

all_tophits <- qtls %>%
  select(gv, gene, `bonf.p`, type) %>%
  group_by(gene) %>%
  slice_min(`bonf.p`, n=1, with_ties=FALSE)
all_tophits$q <- qvalue(sighits$bonf.p)$qvalues

sighits <- all_tophits %>%
  filter(q<=0.05)

sig.hits <- pull(sighits, gv)

data.stronger = mash_set_data(beta.hat[sig.hits,], se.hat[sig.hits,], V=Vhat, df=df.hat[sig.hits,])

the condition has length > 1 and only the first element will be used

U.pca.stronger = cov_pca(data.stronger, 5)
corrplot(U.pca.stronger$PCA_1, is.corr=F)

corrplot(U.pca.stronger$PCA_2, is.corr=F)

corrplot(U.pca.stronger$PCA_3, is.corr=F)

corrplot(U.pca.stronger$PCA_4, is.corr=F)

corrplot(U.pca.stronger$PCA_5, is.corr=F)

corrplot(U.pca.stronger$tPCA, is.corr=F)

Feed this into the ED algorithm

U.ed.stronger = cov_ed(data.stronger, U.pca.stronger)

corrplot(U.pca.stronger$ED_PCA_1, is.corr=F)
corrplot(U.pca.stronger$ED_PCA_2, is.corr=F)
corrplot(U.pca.stronger$ED_PCA_3, is.corr=F)
corrplot(U.pca.stronger$ED_PCA_4, is.corr=F)
corrplot(U.pca.stronger$ED_PCA_5, is.corr=F)
corrplot(U.pca.stronger$ED_tPCA, is.corr=F)

We are maybe starting to get somewhere, but it's still not quite where we'd like it to be. Also worth noting that this is pretty similar to the original one.

corrplot(U.flash$tFLASH_default, is.corr=F)
corrplot(U.flash$FLASH_default_1, is.corr=F)

corrplot(U.flash$FLASH_default_2, is.corr=F)

corrplot(U.flash$FLASH_default_3, is.corr=F)

corrplot(U.flash$FLASH_default_4, is.corr=F)

corrplot(U.flash$FLASH_default_5, is.corr=F)

corrplot(U.flash$FLASH_default_6, is.corr=F)

corrplot(U.flash$FLASH_default_7, is.corr=F)

corrplot(U.flash$FLASH_default_8, is.corr=F)

corrplot(U.flash$FLASH_default_9, is.corr=F)

Error in corrplot(U.flash$FLASH_default_9, is.corr = F) : 
  Need a matrix or data frame!